Rotationally-invariant slave-boson formalism 
and momentum dependence of the quasiparticle weight 
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We generalize the rotationally-invariant formulation of the slave-boson formalism to multiorbital 
models, with arbitrary interactions, crystal fields, and multiplet structure. This allows for the study 
of multiplet effects on the nature of low-energy quasiparticles. Non-diagonal components of the 
matrix of quasiparticle weights can be calculated within this framework. When combined with 
cluster extensions of dynamical mean-field theory, this method allows us to address the effects of 
spatial correlations, such as the generation of the superexchange and the momentum dependence 
of the quasiparticle weight. We illustrate the method on a two-band Ifubbard model, a Hubbard 
model made of two coupled layers, and a two-dimensional single-band Hubbard model (within a 
two-site cellular dynamical mean-field approximation). 
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I. INTRODUCTION AND MOTIVATIONS 
A. General motivations 

The method of introducing auxiliary bosons in order 
to facilitate the description of interacting fermionic sys- 
tems is an important technique in theoretical many-body 
physics. In this regard, the so-called slave boson (SB) 
approach is a very useful tool in dealing with models of 
strongly correlated electrons. Slave boson mean-field the- 
ory (SBMFT), i.e., at the saddle-point level, is the sim- 
plest possible realization of a Landau Fermi liquid (for 
a review, see e.g. Ref.l). Within SBMFT, a simplified 
description of the low-energy quasiparticles is obtained, 
while high-energy (incoherent) excitations are associated 
with fluctuations around the saddle point. In particu- 
lar, two essential features are captured by SBMFT: (i) 
the Fermi surface (FS) of the interacting system (sat- 
isfying Luttinger's theorem) is determined by the zero- 
frequency self-energy, which is in turn determined by the 
Lagrange multipliers associated with the constraints and 
(ii) the quasiparticle (QP) weight Z is determined by the 
saddle-point values of the slave bosons. Hence SBMFT is 
a well-tailored technique when attempting to understand 
the low-energy physics emerging from more sophisticated 
theoretical tools, such as for example dynamical mean- 
field theory (DMFT), which deals with the full frequency 
dependence of the self-energy. 

In this paper, we are concerned with the construction 
of a slave-boson formalism which is able to deal with the 
two following problems. 

1) In multi-orbital models, handle an arbitrary form 
of the interaction hamiltonian, not restricted to density- 
density terms, and possibly including interorbital hop- 
pings or hybridizations. We aim in particular at describ- 
ing correctly the multiplets (eigenstates of the atomic 
hamiltonian), but we also want to be able to work in 



an arbitrary basis set, not necessarily that of the atomic 
multiplets (and of course, to obtain identical results, in- 
dependent of the choice of basis). 

2) Describe situations in which the QP weight is not 
uniform along the Fermi surface, but instead varies as a 
function of the momentum, i.e., Z=Z(k). 

There are clear physical motivations for addressing 
each of these issues. The first one is encountered when- 
ever one wants to deal with a specific correlated material 
in a realistic setting (see e.g. Ref. 2). Usually, more than 
one band is relevant to the physics (e.g. a t2g triplet 
or Cg doublet for transition metal-oxides, or the full 7- 
fold set of /-orbitals in rare-earth, actinides and their 
compounds). The second issue is an outstanding one 
in connection with cuprate superconductors. In those 
materials, a strong differentiation in momentum space 
is observed in the "normal" (i.e., non-superconducting) 
state, especially in the underdoped regime (for a review, 
see e.g. Ref. 3). For momenta close to the nodal re- 
gions, i.e., close to the regions where the superconduct- 
ing gap vanishes, reasonably long-lived QPs are found. In 
contrast, in the antinodal directions, the angle-resolved 
photoemission (ARPES) spectra reveal only a broad line- 
shape with no well-defined QPs. The nature of the in- 
cipient normal state in the underdoped regime (i.e., the 
state achieved by suppressing the intervening supercon- 
ductivity) has been a subject of debate. One possibil- 
ity is that QPs would eventually emerge at low-enough 
temperature in the antinodal region as well, but with a 
much smaller QP weight Zan^Zjsi. Another possibility 
is that coherent QPs simply do not emerge in the antin- 
odal region. Anyhow, there is evidence from ARPES 
and other experiments'^ that the QP weight (whenever 
it can be defined) has significant variation along the FS 
and is larger at the nodes. Since the QP weight sets the 
scale for the coherence temperature below which long- 
lived QPs form, a smaller Z means a smaller coher- 
ence temperature. Hence, if the temperature is higher 
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than the coherence scale associated with momenta close 
to the nodes, and larger than the one associated with 
the antinodes, QPs will be visible only in the nodal re- 
gions. At this temperature, the FS will thus appear as 
being formed of "Fermi arcs" , as indeed observed experi- 
mentally'''. Important differences between the nodal and 
antinodal region in the superconducting state have also 
been unraveled by recent experiments, in particular from 
Raman scattering which revealed two different energy 
scales with different doping dependence, associated with 
each of these regions^. Momentum-space differentiation 
of QP properties is therefore a key feature of cuprate su- 
perconductors, but it is also an issue which is particularly 
difficult to handle theoretically. 

As we now explain, these two issues are actually closely 
related one to the other. In a general multi-orbital model, 
the self-energy is a matrix E^^ (a, (3 are orbital indices). 
Except when a particular symmetry dictates otherwise, 
this matrix has in general off-diagonal (interorbital) com- 
ponents and these off-diagonal components may have 
a non-zero linear term in the low-frequency expansion, 
hence yielding non-diagonal components of the matrix of 
QP weights defined as: 



— : 

du 



(1) 



On the other hand, a momentum dependent QP weight 
Z(k) means that, in real-space, Zij=Z{'R.i — Ilj) depends 
on the separation between lattice sites (a momentum- 
independent Z means that Zij=Z 6ij is purely local). 
Hence, in both cases, one has to handle a QP weight 
which is a matrix in either the orbital or the site in- 
dices. The connection becomes very direct in the frame- 
work of cluster extensions of DMFT (for reviews, see e.g. 
Refs [6-9]). There, a lattice problem is mapped onto a 
finite-size cluster which is self-consistently coupled to an 
environment. This finite-size cluster can be viewed as a 
multi-orbital (or molecular) quantum impurity problem, 
in which each site plays the role of an atomic orbital. Re- 
cently, numerical solutions of various forms of cluster ex- 
tensions to the DMFT equations for the two-dimensional 
Hubbard model have clearly revealed the phenomenon 
of momentum-space differentiation^*^"^^ . Developing low- 
energy analytical tools to interpret, understand, and gen- 
eralize the results of these calculations is clearly an im- 
portant and timely issue. The slave-boson methods de- 
veloped in the present work are a step in this direction. 

Obviously, the existence of off-diagonal components of 
the .Z-matrix is a basis-set dependent issue. A proper 
choice of orbital basis can be made, which diagonalizes 
this matrix. In certain cases, this basis is dictated by 
symmetry considerations, while in the absence of sym- 
metries, the basis set in which Z is diagonal cannot be 
guessed a priori. For instance, in a two-site cluster or 
two-orbital model in which the two sites play equivalent 
roles, even and odd combinations diagonalize not only 
the Z-matrix, but in fact the self-energy matrix itself 
for all frequencies (see Sec. III). In such cases, it may 



be favorable to work in this orbital basis set, and deal 
only with diagonal QP weights. However, performing the 
rotation into this orbital basis set will in general trans- 
form the interacting hamiltonian into a more complicated 
form. For example, starting from a density-density inter- 
action, it may induce interaction terms which are not of 
the density-density type (i.e., involve exchange, pair hop- 
ping, etc.). For these reasons, it is essential to consider 
slave-boson formalisms which can handle both arbitrary 
interaction terms, and non-diagonal components of the 
QP weight matrix: these two issues are indeed connected. 
The formalism presented in this article builds on earlier 
ideas of Li, Wolfle, Hirschfeld and Frcsard^'^'^^ (see Ap- 
pendix. A), in which the SB formalism is formulated in a 
fully rotationally-invariant manner (see also Refs [15,16] 
in the framework of the Gutzwiller approximation), so 
that the orbital basis set needs not be specified from 
the beginning, and the final results are guaranteed to 
be equivalent irrespectively of the chosen basis set. 



B. Some notations 

In this paper, we shall consider multi-orbital models of 
correlated electrons with hamiltonians of the form: 



H 



with 



kin 



^^loc[«] 



k aP 



'k/3 



(2) 
(3) 



In these expressions, a, (3 label electronic species and run 
from 1 to M (i.e., M is twice the number of atomic or- 
bitals in the context of a multi-orbital model of electrons 
with spin: a=(rn, a),<j=\, [). The k vector runs over the 
Brillouin zone of the lattice, whose sites are labeled by 
i (in the context of clustcr-DMFT, i will label clusters 
and runs over the superlattice sites, thus k runs over the 
reduced Brillouin zone of the superlattice, see Sec. HI C). 
The first term in (2) is the kinetic energy: £0,^3 (k) is the 
Fourier transform of the (possibly off-diagonal) hoppings 
and does not contain any local terms (i.e., £a^(k)=0). 
H\oc contains both the one-body local terms and the in- 
teractions, assumed to be local. A general form for H\oc 



isS^: 



a/3 



H\oc — ^ e'^^ndid 



r 4. 

ap-fSdld^dgd^ 



(4) 



Fock states form a convenient basis set of the local 
Hilbert space on each site. They are specified by se- 
quences n=(ni, . . . , um), with na=0,l (we consider a 
single site and drop the site index): 



H = {dl) •••(<) jvac) 



(5) 



In the following, {|^)} will denote an arbitrary basis set 
of the local Hilbert space, specified by its components on 
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the Fock states: 

\A)=J2{n\A)\n) , (6) 

n 

while |r) wiU denote the eigenstates of the local hamil- 
tonian, i.e., the 'atomic' multiplets such that: 

Hioclr) = Er |r) . (7) 



C. Slave bosons for density-density interactions: a 
reminder 

When the orbital densities t?,q are good qnantnm num- 
bers for the local hamiltonian iJioc, i.e., when the eigen- 
states of the latter are labeled by Ua, a very simple slave 
boson formalism can be constructed which is a direct 
multi-orbital generalization of the 4-boson scheme intro- 
duced by Kotliar and Ruckenstein^^'-'^*. While this is 
standard material, we feel appropriate to briefly remind 
the reader of how this scheme works, in order to con- 
sider generalizations later on. We thus specialize in this 
subsection to a local hamiltonian of the form: 

SO that the multiplets are the Fock states |n) themselves, 
with eigenenergies: 

En = ^ Safla + ^ C/a/3 «a . (9) 

To each Fock state, one associates a boson creation op- 
erator Furthermore, auxiliary fermions are intro- 
duced which correspond to quasiparticle degrees of free- 
dom. The (local) enlarged Hilbert space thus consists 
of states which are built from tensor products of a QP 
Fock state, times an arbitrary number of bosons. In con- 
trast, the physical Hilbert space is generated by the basis 
set consisting of the 2^'^ states which contain exactly one 
boson, and in which this boson matches the QP Fock 
state. Thus, the states representing the original physical 
states (5) in the enlarged Hilbert space, in a one-to-one 
manner, are the following ("physical") states: 

|n) = 4|vac)® |n)/ . (10) 

The underlining in |n ) allows one to distinguish between 
the original Fock state of the physical electrons |n), and 
its representative state in the enlarged Hilbert space. In 
this expression, \n)f stands for the QP Fock state: 

\n)f^[fl) ■■■(fl) |vac) . (11) 

It is easily checked that a simple set of constraints 
uniquely specifies the physical states among all the states 



of the enlarged Hilbert space, namely: 

n 

Ena4<^„ = , Va . (13) 

n 

The first constraint imposes that only states with a single 
boson are retained, while the second one insures that the 
fermionic (QP) and bosonic contents match. Obviously, 
the saddle-point values of the slave bosons will have a 
simple interpretation, \(f>„\'^ being the probability associ- 
ated with the Fock space configuration n. 
The operator: 

di=Y.{n\fl\m}cl^i<l^mfl (14) 

nm 

is a faithful representation of the physical electron cre- 
ation operator on the representatives (10), namely: 

dUn)=J2i^'\dl\n)\n') , (15) 

n' 

in which, in fact, the r.h.s is either zero (if na= 1) or 
composed of just a single state (with n'^= 1 and otherwise 
nf3=n'ij for (5^a). This expression of the physical electron 
operators is not unique however: obviously, one can for 
example multiply this with any operator acting as the 
identity on the physical states. This is true as long as the 
constraint is treated exactly. When treated in the mean- 
field approximation however, (i.e., at saddle point), these 
equivalent expressions will not lead to the same results. 
In fact, (14) suffers from a serious drawback namely it 
does not yield the exact non-interacting ({7a/3=0) limit 
at saddle-point. Instead, the expression: 

dl = Y.(^\fl\m} [A„]-^t [1 - A„]-^ /t (16) 

nm 

with A„[,^] = E"«4'^„ , (17) 

n 

turns out to satisfy this requirement, while having ex- 
actly the same action as (14) when acting on physical 
states. This choice of normalization is actually very natu- 
ral given the probabilistic interpretation of jf/inp: the ex- 
pression [Aq]~-^/^0Jj (resp. 0„ [1 — Aq]~-^/^) is actually a 
probability amplitude, normalized over the restricted set 
of physical states such that na=l (resp. na=0). Hence 
the combination of boson fields in (16) is a transition 
probability between the state m with ma=0 and the state 
n with na=l. 

Anyhow, whether the simplest expression (14) or the 
normalized expression (16) is chosen for the physical op- 
erator, the relation between the physical and QP single- 
particle operators is of the form: 

da = fa[<l>]fa (18) 
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It is important to note that the orbital index carried by 
the physical operator is identical to that of the QP oper- 
ator. An immediate consequence is that the self-energy 
at the saddle-point level is a diagonal matrix in orbital 
space T,af3=Saf3'^a, whlch reads: 

= S,(0)+w (^"^) ' ^^^^ 
with Zc = I rap (20) 

Sa(0) = Aa/|ra|2-4 . (21) 

In these expressions, is evaluated at saddle-point level, 
and Xa is the saddle-point value of the Lagrange multi- 
pliers enforcing the constraint (13). 

The expression (20) of the QP weight is an immediate 
consequence of (18): at saddle-point level, becomes 
a c- number and (18) implies that the physical electron 
carries a spectral weight \ra\^. Hence, in order to describe 
within SBMFT situations in which the QP weight is a 
non-diagonal matrix, one must disentangle the orbital 
indices carried by the physical electron and those carried 
by the QP degrees of freedom. These operators will then 
be related by a non-diagonal matrix: 

This is precisely what the formalism exposed in this ar- 
ticle achieves. The physical significance of such a non- 
diagonal relation is that creating a physical electron in 
a given orbital may induce the creation of QPs in any 
other orbital. Thinking of orbital as real-space indices 
(within e.g. cluster-DMFT), this means that the cre- 
ation of a physical electron on a given site induces QPs 
on other sites in a non-local manner, corresponding to a 
momentum-dependent Z{\s.). 



D. Difficulties with naive generalizations to the 
multi-orbital case 

Let us come back to the general multi-orbital inter- 
action (4). In order to motivate the fully rotationally- 
invariant formalism exposed in the next section, let us 
point out some difficulties arising when attempting to 
generalize the simple SB formalism of the previous sec- 
tion. 

The central difference between the general interaction 
(4) and the density-density form (8) is that the atomic 
multiplets |r) are no longer Fock states. Thus, it would 
seem natural to associate a slave boson t^r to each of the 
atomic multiplets. Indeed, Biinemann et al.^^ (see also 
[15]) have proposed generalized Gutzwiller wave func- 
tions in which a variational parameter (a.k.a a probabil- 
ity |?!'rP) is associated with each atomic multiplet (see 
also Ref. [20] and the recent work of Dai et aJ.^^ in the 
SB context). A slave-boson formulation requires a clear 
identification of the physical states within the enlarged 



Hilbert space. A natural idea is to define those in one- 
to-one correspondence with the atomic multiplets, as: 

IE) = 0V|vac) ® ^(n|r) \n)f . (23) 

n 

The local part of the hamiltonian has a simple repre- 
sentation on these physical states i?ioc=X]r -^r'/'r'/'r- 
However, a major difficulty is that there is no simple 
constraint implementing the restriction to these physi- 
cal states, and such that it is quadratic in the fermionic 
(QP) degrees of freedom (which is essential in order to 
yield a manageable saddle point). In particular, it is eas- 
ily checked that the apparently natural constraint^"'^: 

4/a = E<rK|r)</'Ur (24) 
r 

is actually not satisfied by the states (23) as an operator 
identity^^. Further difficulties also arise when attempt- 
ing to derive an expression for the physical creation op- 
erators. These difficulties stem from the fact that two 
atomic multiplets having particle numbers differing by 
one unit cannot in general be related by the action of a 
single-fermion creation. 

One might also think of defining the physical states in 
correspondence to the Fock states, as: 

\n) = \n}f(8)Y^{T\n)cP^^\v<ic) (25) 
r 

which do satisfy the following quadratic constraint: 

4/a = E(rK|r')4'/'r' • (26) 

pp/ 

However, another difficulty then arises. Namely, it is 
not possible to write the local interaction hamiltonian 
purely in terms of bosonic degrees of freedom, which is 
the whole purpose of SB representations. In particular, 
the obvious expression ^?ioc=X^p EY(j>\(j>Y which had the 
correct action on states (23) no longer works for states 
(25) since it leaves unchanged the fermionic content of 
them. 

After some thinking, one actually realizes that these 
naive generalizations are all faced with the same problem, 
namely that they do not embody the crucial conceptual 
distinction between physical and QP degrees of freedom. 
Both (23) and (25) assume a priori a definite relation be- 
tween the physical and QP content of a state. The key 
to a successful SB formalism is therefore to disentangle 
physical and quasiparticle degrees of freedom, and letting 
the variational principle at saddle point decide which re- 
lationship actually exists between those. 

We shall see however in Sec. Ill A that, provided the lo- 
cal hamiltonian has enough symmetries, the rotationally- 
invariant formalism of the present article does correspond 
to assigning at saddle point a probability to each atomic 
configuration (multiplet) [F), hence establishing contact 
with the previous works of Refs. [19,21]. Yet for less sym- 
metric hamiltonians, the general formalism of the present 
article is requested. 
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II. ROTATIONALLY-INVARIANT 
SLAVE-BOSON FORMALISM 

A. Physical Hilbert space and constraints 

In order to construct a SB formalism in which physi- 
cal and QP states are disentangled, we shall associate a 
slave boson (j)-p^ to each pair of atomic multiplet |r) and 
QP Fock state |n)/. More generally, we can work in an 
arbitrary basis set \A) of the local Hilbert space, not nec- 
essarily that of the atomic multiplets, and consider slave 
bosons <?^^„. As we shall see, the formalism introduced 
in this article is such that two different choices of basis 
sets are related by a unitary transformation and there- 
fore lead to identical results. In particular, one could also 
choose the physical Fock states \m)d as the basis set A, 
and work with slave bosons (f)^^ which form the com- 
ponents of a density matrix connecting the physical and 
QP spaces. It is crucial however to keep in mind that the 
first index (A) refers to physical-electron states, while the 
second one (n) refers to quasiparticles. 

A priori, a slave boson can be introduced for any 
pair (A, n). However, in this paper, we shall restrict our- 
selves to phases which do not display an off-diagonal su- 
perconducting long-range order, and hence one can re- 
strict the (py^n's to pairs of states which have the same 
total particle number on a given site (the local hamilto- 
nian i^ioc commutes with J^a^a^^'^)- The formalism is 
easily extended to superconducting states-'^'*'-'^^'^^ by lift- 
ing this assumption and modifying appropriately the ex- 
pressions derived in this section. In the following, we 
consider basis states A which are eigenstates of the local 
particle number (denoted by Na), and hence a is 
introduced provided na=NA- 

The representation of such a basis state in the enlarged 
Hilbert space is defined as: 



1^) 



1 



(27) 



In this expression. Da denotes the dimension of the sub- 
space of the Hilbert space with particle number identical 
to that of A, i.e., Da= D{Na) = {I^J. This insures a 
proper normalization of the state. As before, the "under- 
line" in I A) allows to distinguish this state, which lives in 
the tensor product Hilbert space of QP and boson states, 
from the physical electron state \A). 

Having decided on the physical states, we need to iden- 
tify a set of constraints which select these physical states 
out of the enlarged Hilbert space in a necessary and suf- 
ficient manner. It turns out that the following (M^ -|- 1) 
constraints achieve this goal: 



E^L^^n = 1 (28) 
An 

T.T.^An'^An{n\fiUn') = fife' ,Va .(29) 



A nn' 



physical states arc single-boson states. It is easy to check 
that the physical states satisfy the second set of con- 
straints (29), but a little more subtle to actually prove 
that this set of constraints is sufficient to uniquely select 
the physical states (27) in the enlarged Hilbert space. 
The detailed proof is given in Appendix B. Let us em- 
phasize that the order of primed and unprimed indices 
in (29) is of central importance. 



B. Representation of the physical electron 
operators 

We now turn to the representation of the physical elec- 
tron creation operator on the representatives (27) of the 
physical states in the enlarged Hilbert space. We need to 
find an operator which acts on these representatives ex- 
actly as acts on the physical basis \A). Namely, given 
the matrix elements {A\dl\B) such that 



di\B)=Y,{A\di\B)\A) , 



(30) 



we want to find an operator ^ (in terms of the boson 
and QP operators) such that 

A\B) = Y.^^\d^c.\B)\A) . (31) 



1. Proximate expression 

As in the case of the density-density interactions dis- 
cussed above (Sec. IC), the answer is not imique. We 
first construct the generalization of expression (14) to 
the present formalism (i.e., ignore at first the question of 
the proper operators to be inserted in order to recover the 
correct non-interacting limit). The following expression 
is shown in Appendix C to satisfy (31): 



E 



l3,AB,nm, 



{A\dl\B){n\fl\m) 
^Na{M-Nb) 



t 



^An'+'Bm f 13 



(32) 



We note that Na=Nb -|- 1 in this expression can take the 
values 1, . . . , M. 

Hence, we see that within this formalism, the physical 
and QP operators arc indeed related by a non-diagonal 
transformation (22): 



R[(l)]af3 fp 



(33) 



with the _R-matrix corresponding to (32) given by (i?^^ 
denotes the complex conjugate of Rafj)'- 



{A\di\B){n\fl\m) . 



E fwiM — att 

AB,nm V^a{M-Nb) 



'^An'f'Bm 



(34) 



The first constraint is obvious and requires that the 



The action of (32) on physical states, and the proof that 
it satisfies (31) are detailed in Appendix C. 
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2. Improved expression 

The simple expression (32), although having the cor- 
rect action on the physical states, suffers from the same 
drawback than (14) in the case of density-density inter- 
actions. Namely, at saddle-point level (i.e., with the con- 
straint satisfied on average instead of exactly), the non- 
interacting limit is not appropriately recovered. Thus, 
one needs to generalize the improved expression (16) to 
the present rotationally-invariant formalism. However, 
care must be taken to do so in a way which respects 
gauge invariance (i.e., the possibility of making an ar- 
bitrary unitary rotation on the QP orbital indices, see 
Sec. II C). 

We consider the following operators, bilinear in the 
bosonic fields: 



Anm 



(35) 



(36) 



Anm 



which can be interpreted as particle- and hole- like QP 
density matrices (note that when the constraint is satis- 
fied exactly: A^'^^=Jc/3— A^). We then choose to modify 
the ii-matrix in the following manner (see Appendix C): 

Rmip- E {A\dl\B){n\f\\m)<l>\^<l>j,^M^^, (37) 

AB,nm,j 



with M 



7/3 = \ 7 



/3)(38) 



We chose to let the QP density matrices enter the M- 
matrix in a symmetrized way in order to respect equiva- 
lent treatment of particles and holes. Expression (37) can 
be shown to be gauge- invariant, and turns out to yield 
the correct non- interacting limit at saddle point. How- 
ever, although it yields a saddle point satisfying all the 
appropriate physical requirements, it is not fully justified 
as an operator identity. 



C. Gauge invariance 

As usual in formalisms using slave particles, a gauge 
symmetry is present which allows one to freely rotate the 
QP orbital indices, independently on each lattice site. 
Physical observables are of course gauge- invariant. Let 
us consider an arbitrary SU (M) rotation of the QP op- 
erators: 



4 = E^«/3/i 



(39) 



/3 



This rotation induces a corresponding unitary transfor- 
mation U{U) of the QP Fock states |n)f. This unitary 

transformation is characterized by the fact that the ex- 
pectation value of fa in its Fock basis is an invariant 



tensor: it is the same in every basis. Therefore (summa- 
tion over repeated indices is implicit everywhere in the 
following): 

(n|4|m) = Upp>U{Ur„„,{n'\fl,\m')U{UUn' (40) 

{n\flU\m) = Uaa'U*pf3'^*{U)nn' {n' \fl f ^,\m') U{UU^ 

(41) 

(the second expression can be deduced from the first us- 
ing closure relations) . We can now check that if the slave 
bosons transforms like 



(j)An =U{U)nn'(l>An' 



(42) 



then the constraints and the expressions of the phys- 
ical electron operator (either (32) or (37)) are gauge- 
invariant. Namely, the ii-matrix obeys the following 
transformation law: 

R[4>]ccl3 = R[$\a0' C//3/3' (43) 

and therefore the physical electron operator is invariant: 
da = R[^Ufi3 = R[<PUU (44) 

D. Change of physical and quasiparticle basis sets 

It is clear that the basis \A) of the local Hilbert space 
(i.e., the physical basis states) can be chosen arbitrarily 
in this formalism. Indeed, making a basis change from 
\A) to \A), all the expressions above keep an identical 
form provided the bosons corresponding to the new basis 
are defined as: 



4„ = E<^i^)'^L 



(45) 



As mentioned above, it is often convenient to use the 
eigenstates jF) of iJioc as a basis set. 

Changing the basis states associated with quasiparti- 
cles is a somewhat trickier issue. Up to now, we have 
worked with Fock states |n)/. A different basis set \Q) f 
can be used, provided however the unitary matrix {Q\n) 
is real, i.e., {Q\n) = {n\Q). Indeed, the matrix element 
{Q\n) appears in the transformation of the physical states 
and of the constraint, while {n\Q) appears in the trans- 
formation of the physical electron operator. When this 
matrix elements are real, new bosons can be defined in 
the transformed QP basis according to: 



^'aq 



E(QI")'^An > ((0|n) = (n|Q)) 



(46) 



In particular, when the local hamiltonian is a real sym- 
metric matrix, the same linear combinations of Fock 
states which define the atomic multiplets |r) can be used 
for QPs, and bosons (prr' can be considered. This is 
sometimes a useful way of interpreting the formalism and 
the results at saddle point (see Sec. Ill A). 
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E. Expression of the hamiltonian, free energy and 
Green's function 

In this section, wc derive the expression of the hamil- 
tonian in terms of the slave boson and QP fermionic vari- 
ables. We then construct the free-energy functional to be 
minimized within a mean- field treatment, and express the 
Green's function and self-energy at saddle point. 

We recall that the full hamiltonian (2) reads, in terms 
of the physical electron variables: i?=-ffkin + X^i-^iocW 
with -ffkin=ZlkZla/3^"/3('*) the intersite kinetic 
energy and H\oc the local part of the hamiltonian on a 
given site i, with general form (4). 

It is easily checked that the following bosonic operator 
is a faithful representation of i?ioc on the representatives 
of the physical states in the enlarged Hilbert space: 



'Bn 



(47) 



AB 



If the basis |r) of atomic multiplets is used, this simplifies 
down to: 



(48) 



Using the bosonic i?-operators relating the physical elec- 
tron to the QP operators, yields the following expression 

of the kinetic energy: 



H 



kin 



[R^U'ea'0'm'0fLfy.is (49) 



A mean-field theory is obtained by condensing the slave 
bosons into c-numbers {(pAn) = 'PAn- The constraints 
are implemented by introducing Lagrange multipliers: Ao 
associated with (28) and A^q' = [AJaa/ associated with 
(29). The saddle point is obtained by extremalizing, over 
the ^pAnS and the Lagrange multipliers, the following 
free-energy functional: 



n[{(^A«};A,Ao] 



(50) 



-/3(Rt(<^)£(k)R(<^)+A) 



Ao 



+ XI '^*An' { ^uti'Sab Ao -I- 5„„'{A\Hioc\B) 

ABnn' 

-<5ab E A«/3(n|/J,//3|n') } (^i3„ 

a/3 



(51) 



The saddle-point equations, as well as technical aspects 
of their numerical solution, are detailed in Appendix D. 

Finally, we derive the expressions of the Green's func- 
tions G, the self-energy S and the QP weight Z at sad- 
dle point. For the QPs, the one-particle Green's func- 
tion Gf^a^iK T - T')=-{flcxi'^)hl3{T')) reads (in matrix 
form) : 



Rt(v5)£(k)R(<^) - A (52) 



and hence the physical electron Green's function reads 
(we drop the ip dependence for convenience): 

G^i(k,a;) = [Rt]-iG7iR-i 

= a;(RRl')-i - [RI'J-IAR-i - £(k) , (53) 

while the non-interacting Green's function is (including 
the one-body term present in i?ioc): 



G^o^(k,a;) =a;l-e° -£(k) 



(54) 



The physical self-energy is thus: 

S,(c.) ^ G-J G-/ 

= Lu{l- [RRt]-i) + [Rt]-iAR-i - £° . (55) 

So that the matrix of QP weights is obtained in terms of 
the .R-matrix at saddle point as: 



Z = RR 



(56) 



This generalizes (20) to non-diagonal cases. It is easily 
checked that these expressions of the physical quantities 
Gd, and Z are indeed gauge- invariant. 



III. ILLUSTRATIVE RESULTS 

In the following, we apply the above formalism to three 
different model problems in strongly correlated physics. 
First, we consider two popular models, namely the two- 
band Hubbard model on a three-dimensional (3D) cubic 
lattice, and a "bi-laycr" model, coupling two Hubbard 3D 
cubic lattices . Finally a two-site cluster (cluster-DMFT) 
approximation to the single-band Hubbard model on 
a two-dimensional (2D) square lattice is investigated. 
Hence these models have in common that they all involve 
two coupled orbitals (associated, in the cluster-DMFT 
(CDMFT) framework, to the dimer made of two lattice 
sites). The present formalism is of course not restricted 
to two-orbital problems, however such models provide the 
simplest examples where the power of the method may 
be demonstrated. 



A. Two-band Hubbsird model 

The Hubbard model involving two correlated bands, 
without further onsite hybridization or crystal-field split- 
ting, serves as one of the standard problems in condensed 
matter theory. In contrast to the traditional single-band 
model, the formal interaction term in Eq. (4) now gener- 
ates in the most general fully SU (2) symmetric case four 
energy parameters, i.e., the intraorbital Hubbard U and 
the interorbital Hubbard U' as well as the two exchange 
couplings J,Jc- Thus the present atomic hamiltonian 
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TABLE I: Eigenstates |r) of the SU{2) rotationally-invariant two-band Hubbard model. Spin values and energies are given for 
the eigenstates. The last column shows the slave bosons for the description of the eigenstates in the SBMFT formalism. 



No. 




Eigonstate |r} 


Sr 




Er 




4>rn 


1 




|00,00) 













•^1,100, 00) 


2 
3 
4 
5 




1 TO, 00) 
|0i,00> 

|oo, T 0) 

|00,0i) 


1 

2 
1 
2 
1 
2 
1 
2 


1 

2 

1 

"2 
1 
2 
1 

"2 










'i^2,|T0,00> 

'A3,|o;,oo> 

"^4,100,10) 
4'5,\00,Ol) 


6 
7 
8 
9 
10 
ii 


1 T 0, T 0) 

;75(l T0,0i) + |0i,T0)) 
|0i,0i> 

Tfd TO,oi)-loi,TO)) 
^(1 Ti,oo)-|oo,n)) 
^(1 Ti,oo) + |oo,n)) 


1 
1 
1 




U 


1 



-1 




U 


U'-J 
U'-J 
U'-J 

U' + J 
U-Jc 

U + Jc 




'?^6,1T0,T0) 
('/'r.iTO.oi)' "^T.ioi.TO)) 
't's,\oi,oi) 

(</'9,|T0,0i>> <A9,|0J.,T0>) 

('^io,|Ti,oo>> "T^io.ioo.n)) 

( ± ± \ 

Vv'ii.in.oo) ' v'lijoo.n)^ 


12 
13 
14 
15 




1 Ti,TO) 

in,oi> 

1 To,n) 

|Oi,Ti) 


1 
2 
1 
2 
1 
2 
1 
2 


1 
2 
1 

"2 
1 
2 
1 

"2 


U + 2U' - 
U + 2U' - 

U + 2U' - 
U + 2U' ~ 


J 
J 

J 
J 


?^l2,|Ti,T0) 
</'l3,|Ti,0J.) 

'^i4,iTo,n> 


16 




1 n,n> 








2U + AU' - 


2J 





reads 

a crcr' 

-j'^nian2a + J'^d\^dl^dl^d2a 

a cr 

+Jc {d\^d\^d2id2-[ + d\-^d\^d\idi-^ . (57) 

The kinetic energy shall contain only intraband terms 
for a basic tight-binding (TB) model for s-bands on a 
3D simple cubic lattice with lattice constant a. Thus the 
corresponding hamiltonian is written as 

^?kin = - i E E E '^IcJjaa > (58) 
<T a=l,2 i^j 

with the eigenvalues 

ea(k)=--fa ^ cos{kfj,a) , (59) 

IJ.=xyz 

where ta denotes the hopping parameter for orbital 
a=l,2. For convenience, we set o=l. The factor 1/3 in 
eq. (58) is to normalize the total bandwidth to Wa=^ta- 
Because of the cubic symmetry, U'=U — 2 J may be used, 
and furthermore we set J=Jc- This model is similar to 
the one considered by Biinemann et al}^ using a gener- 
alized Gutzwiller approximation (see also [15]). Our sim- 
pler TB description exhibits in principle perfect nesting, 
however this issue is not relevant at the present level. The 
3D two-band Hubbard model is studied to make contact 
with the named previous work and in order to establish 



the connection between the Gutzwiller and slave-boson 
points of view. 

When working in the SU (2) rotationally-invariant 
case, the 2^=16 atomic eigenstates |r) of the local hamil- 
tonian (57) serve as the appropriate atomic basis (see Ta- 
ble I) , however also the simpler Fock basis (or any other) 
may be used. Of course, in the Fock basis, a more com- 
plicated energy matrix must be used in the saddle-point 
equations (see Appendix D). It should be clear from 
Tab. I that there are 20 nonzero slave-boson amplitudes 
0P„ for the current problem. The 5*^=0 triplet as well 
as the three singlets are described with two (/)p„, respec- 
tively. In principle, even more (^p„ may be introduced in 
the beginning of the iteration cycle to minimize f7, but at 
convergence those one will come out to be strictly zero. 
Of course, in high-symmetry situations there is still some 
redundancy within the set of the 20 SBs. For instance, 
for equal bandwidth at half filling (see Fig. 1), all the 
one- and three-particle SBs are identical, as well as the 
zero- and four particle SB. Moreover the S*^ = ±1 triplet 
SBs are equal because of the degeneracy. The two SBs 
describing the S'^=0 triplet are also identical, with a mag- 
nitude <i>(^to)n~^{t ±i)nl^- bosons describing 
one specific singlet have the same absolute value, however 
they carry the multiplct phase information, i.e., have plus 
or minus sign. In conclusion, in the orbitally degenerate 
case, the SB amplitudes at saddle point are of the form: 

<^r„ = (n|r)yr , (60) 

in which the matrix element (r|n) is entirely determined 

by -ffioc and ?/r is a (coupling-dependent) amplitude, de- 
pending only on the eigenstate F. This is more clearly in- 
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0.3 



0.2 



0.1 



0, 



1 

'l'l,IOO,00> 


1 

1 

1 


^2, 110,00> 




- — 'I'lO, lll.OO> 




'I'e, iio,io> 




Til, lll,00> 















0.5 1.0 1.5 2.0 

U/(2t) 



FIG. 1: (color online) Inequivalent Slave-boson probabili- 
ties |<^rnP for the two-band Hubbard model at half-filling 
for equal bandwidth and J/C/=0.2. Note that (j>io m oo> ^^^'^ 
"/"ii.lT-l.oo) describe part of the singlet states, hence their over- 
all amplitude is scaled by l/y/2. 




FIG. 2: (color online) Influence of J on the Mott 
transition in the two-band Hubbard model at half fill- 
ing (n=2) for equal bandwidth. From right to left: 
J/U=0, 0.01, 0.02, 0.05, 0.10, 0.20, 0.45. 



terpreted when atomic states are also used as basis states 
for QPs (Sec. II D). Indeed, Eq. (60) means that: 



(61) 



Hence, m this highly symmetric case, the saddle point is 
indeed of the diagonal form considered in Refs. [19,21]. 
Once the symmetry is lowered, more SBs become inequiv- 
alent and this relation does not hold anymore: there are 
off-diagonal components even when the basis of atomic 
states is used for both physical and QP states. In this 
context, the present formalism becomes essential. Dif- 
ferent bandwidths for each orbital, together with a finite 
doping away from half-filling lead for instance to two dif- 
ferent absolute values for the two SBs associated with the 
singlets formed by the two doubly-occupied Fock states 
(as seen at the end of this paragraph in Fig. 6). 

Since no interorbital hybridization is applied in this 
section, the Z-matrix is diagonal. We consider first the 
simple case of equal bandwidths ti=t2=Q-5 (note that 
in all our applications, t sets the unit of energy), thus 
Z\\=Zi2=Z . Figure 2 shows the variation of Z for 
different ratios J jU in the half-filled case (n=2). The 
critical coupling Uc for the Mott transition with J=0 ob- 
tained from this slave-boson calculation is in accordance 
with the result of the analytical formula given by Fresard 
and Kotliar^*. It is seen that an increased J lowers the 
critical U and moreover changes the transition from sec- 
ond to first order. Note that in this regard. Fig. 2 depicts 
Z up to the spinodal boundary, i.e., the true transition 
(following from an energy comparison) is expected to be 
at slightly lower Uc- One can also observe the nonmono- 
tonic character for the evolution of the critical Z at this 
boundary when increasing J jU . We plot in Fig. 2 addi- 
tionally the results when restricting the atomic hamilto- 
nian to density-density terms only, in order to check for 
the importance of the then neglected spin-flip and pair- 
hopping terms. For larger J jU the critical Z from the 



latter description is larger compared to the rotationally- 
invariant one and moreover it is monotonically grow- 
ing. The latter feature strengthens the first-order char- 
acter in the density-density formulation for growing J/J7, 
whereas for rotationally-invariant interactions this char- 
acter is strongly weakened in that regime. Although for 
J/J7=0.45 the jump of Z is quite small, the transition is 
however still first order in the present calculation. Fur- 
thermore, there appears to be a crossover between the 
two approaches concerning the reachable metallic spin- 
odal boundary when increasing J jV . 

At quarter filling (n—\-) a continuous transition is ob- 
tained for all the previous interactions (see Fig. 3) . Com- 
pared to the half-filled case, the density-density approx- 



rotationally-invariant 
density-density 




J/U=0.20" 



FIG. 3: (color online) Influence of J on the Mott transi- 
tion in the two-band Hubbard model at quarter fllling (n=l) 
for equal bandwidth. The first three combined curves for 
the two types of interactions (from left to right) belong to; 
J/(7=0, 0.05, 0.10. Arrows indicate the labelling for the two 
larger J jU ratios. 
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FIG. 4: (color online) Influence of J for fixed U at n=2 (solid 
lines) and n=l (dotted-dashed lines) for equal bandwidth and 
full SU{2) symmetry. The vertical dotted lines mark the limit 
we set for J, respectively. 



imation appears to be less severe for small J/U, but 
leads to some differences compared to the rotationally- 
invariant form for large J/U. Note that for J/U=0A5, 
U'—J in the local hamiltonian (57) becomes negative. 
Thus a corresponding change of the ground state may 
lead to the resulting nonmonotonic behavior for Uc then 
observed in Fig. 3. The critical U for J—Q is smaller than 
at half filling and with increasing J/U the transition is 
shifted to larger Uc (with the above named exception for 
J/U large). Hence J has a rather different influence on 
the degree of correlation for the two fillings. While for 
n=2 the Hund's rule coupling substantially enhances the 
correlations, seen by the decrease in Z, for n=l the oppo- 
site effect may be observed. This is also demonstrated in 
Fig. 4 which displays the influence of J for fixed values of 
U comparing half filling with quarter filling. The strong 
decrease in Z upon increasing J was recently shown to be 
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1 < 1 , 1 < 1 < 



0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 
filling 

FIG. 5: Filling dependence of Z for selected values of J/U 
within the equal-bandwidth two-band model with full SU{2) 
symmetry {U—1.75). 



important for the physical properties of actinides, in par- 
ticular regarding the distinct properties of (5-Plutonium 
and Curium^'^. For each U shown in Fig. 4, the density- 
density limiting value J7/3 was used as an upper bound 
for J. However, for U—2 and U—3 the system shows 
already a first-order transition at half-filling below the 
latter limit. 

The QP residue Z is shown as a function of filling n 
in Fig. 5 for U^l.75 and three ratios J/U. For J=0 it 
is observed that Z{n) exhibits two minima, both located 
at integer filling. The minimum at n—l is deeper, cor- 
responding to a lower value for Uc in the quarter-filled 
case. Because of the filling-dependent effect of J seen in 
Fig. 4, the nonmonotonic character of Z{n) is lifted for 
growing J/U. 

The two-band Hubbard model was already extensively 
studied in the more elaborate DMFT framework in in- 
finite dimensions. Such investigations reveal the same 
qualitative change of the critical U for different integer 
fiUings^'*'^^, of course with some minor quantitative dif- 
ferences. Also the reduction^'''^^"^^ of Uc and the onset 
of a first-order Mott transition^'*'^^'^^^^^ for finite J/U at 
half-filling is in accordance. Concerning the latter effect, 
the trend of weakening the first-order tendency for large 
J/U is also reproduced and there is some discussion^^'^^ 
about the possibility of even changing back to a contin- 
uous Mott transition in that regime. The increasing Uc 
with growing J at quarter filling was also found by Song 
and Zou^^. 

Finally, in Fig. 6 a comparison between the equal- 
bandwidth and the different-bandwidth cases at noninte- 
ger filling n=1.5 is displayed {J/U— 0.2). For W\ = W-2. the 
model does not show a metal-insulator transition because 
of the doping. Also the filling of both bands is iden- 






1 2 3 4 
U/(2t) 

FIG. 6: (color online) Comparison of the two-band model 
for W\=Wn (left) and Wi=2W2 (right) at noninteger filling 
n=1.5 and J/U=0.2. The ratio x is plotted for the singlet 
states coupling the doubly occupied Fock states (see Tab I), 
demonstrating that iprr' is no longer diagonal in this case. 
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tical and constant with increasing U (nsi=ns2=0.375), 
and as stated earlier the SBs are still of the form given 
by eq. (60). However, when breaking the symmetry be- 
tween the two bands by considering different bandwidths, 
the model behaves qualitatively rather differently. The 
individual band fillings are not identical anymore, favor- 
ing the larger-bandwidth band for U=Q. With increas- 
ing U the system manages to drive at least one band 
insulating by transferring charge from the broader into 
the narrower band, until the latter is filled with one elec- 
|-j-Qj^30,3i HQYice Z2 of the narrower band becomes zero at 
an orbital-selective Mott transition (OSMT)3°"35^ r^j^j^ 
asymmetric model has also a more sophisticated SB de- 
scription, since for instance the SBs of the singlets built 
out of the respective doubly-occupied Fock states have 
now different amplitudes. 

B. The Hubbard bilayer 

Next, we consider a model consisting of two single- 
band Hubbard models (two "layers"), coupled by an in- 
terlayer hopping V. This rather popular model has al- 
ready been subject of various studies^®"'^^ . For simplicity 
and in order to make connection to the previous section, 
each layer is described here by a 3D cubic lattice, with 
an onsite repulsion U and an intralayer bandwidth Wa 
{a — 1,2), possibly different for the two layers. Hence 
the local hamiltonian for this problem reads 

Hloc = U na-lUal +VY^ (4o-C^2<T + 4^C^1<t) 

+ ^E'^l^l-'4a'rf2a , (62) 
crcr' 

where the last term describes a possible spin-spin interac- 
tion between the layers. However, for simplicity, we only 




U/(2t) 

FIG. 7: (color online) Half-filled bi- layer with equal band- 
width for V=Q and 1^=0.25. For V=Q the filling per spin 
within the two bands is identical (nsi=ns2=0.5), whereas 
for finite V the the symmetry-adapted bonding/antibonding 
states have different filling, denoted n+,n_. 




>l , I , \ J I , \ , \ , I , L 

0.5 1.0 0.5 1.0 1.5 2.0 

U/(2t) 



FIG. 8: (color online) QP residues Zi and symmetry-adapted 
fillings n+,n- for the half-filled bilayer. Left: W2/W1—O.5 
and V=0. Right: various bandwidth ratios and V—Q.l. In 
the right part, the curves for smaller Z and n are associated 
with the lower-bandwidth band. 



present in this article results with ,7=0. Our choice of 
kinetic energy is equivalent to the one in the last section, 
i.e., given by Eqs. (58,59). 

In the presence of V, an off-diagonal self-energy 
Si2(cij) is generated. Furthermore, away from half-filling 
{ni+n2=2), this self-energy is expected to have a term 
linear in uj at low frequency, and hence ^127^0. We note 
that, when the bandwidths are equal {Wi=W2), the bi- 
layer model can be transformed into a two-orbital model 
by a k-independent rotation to the bonding-antibonding 
(or ) basis. In the latter basis, there is no hybridiza- 
tion but instead a crystal-field splitting (=21/) between 
the two orbitals. The couplings of the two-orbital hamil- 
tonian are given by (in the notation of the previous sec- 
tion, and for J=0): UcS^U^fi^JcS=U/2. When the 
bandwidths are different however, the interlayer hopping 
cannot be eliminated without generating non-local inter- 
dimer interactions. 

Due to the reduced symmetry of the present model 
in comparison to the two-band Hubbard model from the 
previous section, the number of nonzero SBs equals 
now 36 (we use here the Fock basis for \A)). We first con- 
sider the simplest case of a half-filled system {ni—n2—l) 
with equal bandwidths W1—W2 (and J~0). Results for 
the intralayer QP weight and the orbital occupancies of 
the bonding and antibonding bands are given in Fig. 7. 
It is seen that the Mott transition is continuous for V—Q 
but becomes discontinuous in the presence of an inter- 
layer hopping Vy^O. For V=0.25 the spinodal boundary 
of the metallic regime is reached for ?7'^2.055. These re- 
sults are consistent with findings in previous works'^^"'*" 
within the DMFT framework. 

Still focusing on the half-filled case, we display in Fig. 8 
the QP weight as a function of U for different bandwidth 
ratios W2/W1. When V—0, one has two independent 
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).0 0.5 1.0 1.5 2.0 2.5 
U/(2t) 

FIG. 10: (color online) Bilayer at fixed doping (n=1.88) with 
equal bandwidth and 1^=0.25. Full green (gray) lines: fillings 
n+ ,11- , dashed dark lines: QP weights Z+ , Z- . The vertical 
dotted line marks the critical U at half filling. 



Mott transitions in each layer, i.e., an OSMT scenario, 
at which Z vanishes continuously. In the presence of a 
non-zero V, this is replaced by a single discontinuous 
transition for both orbitals. This is consistent with pre- 
vious findings on the OSMT problem'^^. 

We now consider the effect of finite doping away from 
half filling. Fig. 9 displays the diagonal {Zii=Z22) ele- 
ments as well as the now appearing Z12 element of the 
QP weight matrix as a function of doping for U<Uc- Ad- 
ditionally shown are the symmetry-adapted QP weights 

(occupations which follow from diagonalizing 

the Z (A(P)) matrix. Since Z12 is small in this case, the 

are rather similar to .^11=^22 and merge with the 
latter at half-filling. On the other hand, the polarization 
of the (-I-, — ) bands is still increasing. 

As it is seen in Fig. 10 the off-diagonal component Z12 
becomes increasingly important for larger U {>Uc) in the 
doped case. It follows that in this regime the QP weights 

for the bonding/antibonding bands have rather dif- 
ferent magnitude/behavior. Whereas Z- is monotoni- 
cally decreasing, Z^ turns around and grows again (as 



also is the filling of the bonding band). Hence, this model 
is a simple example in which a differentiation between QP 
properties in different regions of the FS occur. Fig. 11 
shows the QP (-1-,— )-bands in the noninteracting and in- 
teracting case {U>Uc), exhibiting strong orbital polariza- 
tion and different band narrowing close to the insulating 
state. For very small doping and large U a transition to 
a new metallic phase is found, which will be discussed in 
detail in a forthcoming publication''^ . 

Finally, we have also investigated a case in which 
the interlayer (interorbital) hopping does not have a lo- 
cal component {V=0), but does have a non-local one 
V=ti2y^0, treated in the band term of the hamiltonian. 
Hence the corresponding energy matrix reads here 



£(k) 



5" J- ] Y: cos(fc,a) 

y fi—xyz 



(63) 



with the choice fii=i22=0.5 and ii2=0.25, as well as 
a=l. In that case, a continuous Mott transition within 
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FIG. 12: (color online) Half-filled bilayer, with equal band- 
widths and V—Q, but with a non-local interlayer hybridization 

tl2. 
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an OSMT scenario can be recovered, with, interestingly, 
a sizeable value of the ofF-diagonal Zyi (Fig. 12). At the 
transition Z\\=Z22=Z\2=Z ^ holds, i.e., the Z-matrix has 
a zero eigenvalue, associated with the (antibonding) in- 
sulating band. Note however that no net orbital polar- 
ization appears with V being purely non-local. 



C. Application to the momentum-dependence of 
the quasiparticle weight within cluster extensions of 
DMFT 

In this section, we finally consider the implications of 
the rotationally invariant SB technique for the Mott tran- 
sition and the momentum-dependence of the QP weight, 
in the framework of cluster extensions of DMFT. 

For simplicity, we consider a CDMFT approach to the 
two-dimensional Hubbard model with nearest-neighbor 
hopping t and a next-nearest neighbor hopping t', based 
on clusters consisting of two sites (dimers), arranged in 
a columnar way on the square lattice (see Fig. 13). The 
"local" hamiltonian on each dimer is formally identical 
to the one introduced in the previous section for the bi- 
layer model, i.e. Eq. (62), with the value V=—t of the 
inter-'orbital' hybridization. The interdimer kinetic en- 
ergy matrix reads (we set again a—X): 



eii(K) = £22(K) = -2t cosKy 
ei2(K) = £*i(K)--ie 



(64) 



2t' (I 



-e'"^^-) cosK,, 



in which K denotes a momentum in the reduced Bril- 
louin zone (BZ) of the superlattice: € [—tt/2,+tt/2], 
Ky G [— TT, -|-7r]. Note again that in SB calculations, the 
intradimer t has to be treated separately from the rest 
of the kinetic energy within Hioc- It is easy to check 
that when putting back —t into the offdiagonal elements 
of the above kinetic-energy matrix, the eigenvalues just 
correspond to the one of a single band: 



e(k) = — 2f(cos kx + cos ky) — At' cos k^ cos ky 



(65) 
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FIG. 14; (color online) Half-filled two-dimensional Hubbard 
model within 2-site CDMFT. Left: QP weights and band 
fillings, right: static dimer self-energy E. 



in the full BZ of the original lattice kx,y G [— tt, -|-7r]. 
In the following, we set t=0.25 and consider successively 
t'=Q and t'= — 0.3i (a value appropriate to hole-doped 
cuprates). Note that, in this article, we do not consider 
a bigger cluster than the 2-site dimer, even in the pres- 
ence of t' . Hence, the cluster self-energy will only con- 
tain Sii and E12 components, i.e., has a spatial range 
limited to the dimer. As a result, no renormalization of 
the effective t' is taken into account. This is of course 
an oversimplification (particularly in view of the demon- 
strated physical importance^^'"'^^ of S13 close to the Mott 
transition) . Larger clusters will be considered within the 
present SBMFT in a further publication. The goal of the 
present (simplified) study is to make a point of princi- 
ple, namely that the SB formalism can indeed produce 
a momentum-dependent ^(k). As the cluster symme- 
try of the problem at hand is identical to the bilayer 
model from the last section, the number of nonzero SBs 
amounts again to 36. Figures 14, 15 and 16 summa- 
rize the main findings, at half- filling and as a function of 
doping, respectively. Let us first concentrate on the case 
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FIG. 13: Square lattice in the 2-site CDMFT picture. 



FIG. 15: (color online) Doped two-dimensional Hubbard 
model [t'^Q) within 2-site CDMFT, for (7=2.195 ((7c~2.197). 
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FIG. 16: (color online) Two-dimensional Hubbard model 
within 2-site CDMFT at fixed doping n=0.94. Left: t-only 
model, right: t-t' model {t' ——0.3t). The vertical lines denote 
the critical U at half fiUing, respectively. 



t'=0. Obviously, the Mott transition at half-filling in this 
case occurs in a manner which is very similar to the bi- 
layer model with a finite interlayer hybridization studied 
in the previous section (Fig. 7): a /irs<-onier transition is 
found. The static part of the self-energy I]ii=S22 equals 
t//2, while S12 (which has no frequency dependence at 
half-filling within SBMFT) has a more complicated nega- 
tive amplitude close to the transition. Note that we only 
discuss the paramagnetic solution, though, of course, the 
system is in principle unstable against antiferromagnetic 
order for any U. Upon doping, a finite value of Z12 is 
generated. The behavior of Z12 is rather similar to the 
case of the bilayer model, except for the change of sign 
(Figs. 15, 16). Hence its amplitude is again significantly 
enhanced for U>Uc, i.e., the values tend to mani- 

festly deviate from each other. This is therefore signalling 
an increasingly nonlocal component of ^(k) as the Mott 
insulating state is approached at strong coupling. Again, 
further studies in the latter regime at small doping will 
be published soon**^. 

Including the effect of a nonzero nearest-neighbor hop- 
ping turns out to lead to significant differences. Al- 
though the first-order character of the transition remains 
stable, the critical U is significantly lower (Fig. 14). The 
static components of the self-energy behave rather simi- 
larly to the i-only case, with some minor quantitative dif- 
ferences. There is a small negative .^12 with a maximum 
amplitude ^--^0.02, remaining nonzero also at the Mott 
transition (~0.01). The main difference in comparison 
to t'=0 is that here, in the doped case, Z12 changes sign 
from negative to positive close to the insulating regime 
for U>Uc+S (with S>0) (see Fig. 16). Thus the degree of 
correlation of the effective (bonding-antibonding) bands 
is inverted. These differences have to be interpreted with 
caution however, since again the hopping range on the 
lattice is larger than our cluster size, and definitive con- 
clusions will have to be drawn from a study involving S13 
as well. 

Nonetheless, keeping with the simplified treatment 



based on a 2-site cluster, we now describe the resulting 
momentum dependence of the QP weight Z(k) for the t-t' 
model. The matrix elements Sn and S12 of the cluster 
(physical) self-energy matrix Sc are obtained from the 
SB amplitudes at saddle point according to (55). The 
self-energy is then periodized on the whole lattice, in the 
form42.43. 

i;iat(k,Lj) = I]ll(w) + i Ei2(w)(cOsfc:r + COS fcy) . (66) 



The interacting FS is defined as follows 

Ai-e(k)-Siat(k,w=0) = . 
For our case, using Eqs. (65, 66), this reads 

1 



(67) 



At-Sn(0) + 



2t 



Si2(0) 



(cos kx + cos ky) 
+4:t' cos kx cos ky = . (68) 



Hence the FS deforms in a nontrivial way in the presence 
of Slat, when including t' in the present 2-site CDMFT 
description. The QP weight ^(k) can be derived from 
Slat according to: 



Z(k) 

which leads here to: 



d 

1 - ^Siat(k,tj) 



k=kF 



Z(k) 



(■^11 ^ ^12) 



12 (cos kx + COS ky) 



Zii — -Zi2{coskx + cos ky) 



(69) 



(70) 

-1 



A contour plot of this function is displayed in Fig. 17. 
Note that it varies only according to (cos kx + cos ky ) . Be- 
cause the interacting FS involves both t and t' , and hence 




FIG. 17: (color online) Interacting Fermi surface (solid lines) 
for the CDMFT treatment of the 2D t-t' Hubbard model with 
t'=—0.3t and U—2.5 at n=0.94 (per site). The color contours 
show the variation of Z(k) (smallest at antinodes). 
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both lattice harmonics (cos + cos ky) and cos k^ cos ky 
(for t'y^O), it cuts through different contour hues of Z(k). 
This results in a QP weight which varies on the FS. Fig- 
ure 17 shows ^(k) for k close to the interacting FS. 
Albeit the momentum variation is quantitatively quite 
small, the key qualitative effect of Z being different on 
different part of the FS is indeed found. It is seen that 
the QPs along the nodal direction, i.e., along (0,0)-(7r, tt), 
have slightly larger Z than the ones in the the antinodal 
direction ((0,0)-(0, tt)). Hence these results are indeed 
in qualitative agreement with ARPES measurements on 
cupratcs. Note that to get nodal points to be more co- 
herent that antinodal ones in this 2-site scheme, Zi2>0 
is actually crucial. 

Our results provide, to our knowledge, the first ex- 
ample of a SB calculation which can address the issue of 
the momentum dependence of the QP weight. We believe 
that the too small variation of Z along the FS found here 
is due to the oversimplified 2-site description in which S13 
is neglected. We intend to consider improvements on this 
issue using the present SBMFT in a forthcoming work. 

Finally, let us make contact with previous work on 
the two-dimensional Hubbard model. Of course, this 
model has been intensively studied with a variety of 
methods such as: quantum Monte-Carlo^^'^^, exact diag- 
onalization^^'^^, path- integral renormalization group***^. 
functional renormalization group**^ and various quan- 
tum cluster methods (dynamical cluster approxima- 
tion^°, cluster extensions of dynamical mean- field the- 
ory'^, and variational cluster perturbation theory^). We 
shall not attempt here a detailed comparison between 
the rotationally-invariant SB method (which anyhow is 
a mean-field technique tailored to address low-energy is- 
sues) with the results of these numerical methods over 
the whole phase diagram (note in particular that we 
have not yet investigated long-range ordered phases, such 
as antiferromagnetism or superconductivity). Rather, 
we would like to point out that some recent numerical 
studies using the above methods^^"^^'^^ have indeed re- 
vealed the emergence of momentum-space differentiation 
in the two-dimensional Hubbard model. We hope that 
the rotationally-invariant SB method will help under- 
stand qualitatively the low-energy physics emerging from 
these results. 



IV. CONCLUSION AND PERSPECTIVES 

In this paper, we extended and generalized the 
rotationally invariant formulation of the slave-boson 
method^^'". 

Our formulation achieves two goals: (i) extending the 
slave-boson method in order to accommodate the most 

general crystal fields, interactions and multiplct struc- 
tures and (ii) the development of a technique which can 
describe QP weights and Fermi liquid parameters which 
vary along the Fermi surface. 

The key aspect of the formalism is to introduce slave 



boson fields which form a matrix with entries labeled by 
a pair of a physical state and a QP state (within an ar- 
bitrary choice of basis set). As a result, a density matrix 
is constructed instead of just a probability amplitude for 
each state. 

While the first objective (i) could also be achieved 

by generalizing appropriately the Gutzwillcr approxima- 
^jQj^i5,i6^ wc find the slave-boson approach to be some- 
what more flexible, in the sense that it is a mean field 
theory which can in principle be improved by computing 
fluctuations around the saddle point. Our application to 
the two band model seems promising. While further work 
is needed to benchmark the accuracy of the rotationally- 
invariant slave-boson method against exact quantum im- 
purity solvers, it is clear that already in the single site 
multiorbital DMFT setting, our method has numerous 
advantages. It obeys the Luttinger theorem even in the 
presence of multiplcts, and can accommodate full atomic 
physics information. Furthermore, the off-diagonal el- 
ements of the matrix of QP weights can be calculated 
within this method, while the standard slave-boson or 
Gutzwiller approximations (using probability amplitudes 
instead of a density matrix) cannot achieve this goal. 

Our technique achieves the second objective (ii) via a 
detour, namely the use of cluster extensions of dynami- 
cal mean-field theory in order to reduce the lattice to a 
multisite (molecular) impurity problem, to which we ap- 
ply our rotationally-invariant slave-boson method as an 
impurity solver. Because the intersite matrix elements 
of the QP weight can be calculated, it leads on the lat- 
tice to a momentum-dependence of the QP residue Z{\s.). 
We successfully demonstrated this point, in the frame- 
work of a 2-site CDMFT study of the single-band 2D 
Hubbard model. We did find that the QP weight at 
the nodes is somewhat larger than at the antinodes, al- 
though the magnitude of this effect is expected to in- 
crease within a more realistic study involving a larger 
cluster (e.g., a square plaquette), which is left for fu- 
ture work. A major challenge is the direct extension of 
our slave-boson approach to the lattice, without resorting 
to the cluster-DMFT detour. In this context, we men- 
tion that other slave-boson techniques, which introduce 
magnetic correlations through the use of link variables to 
decouple the superexchange J term^, can be interpreted 
in terms of a k dependent self-energy. However, within 
such schemes, the derivative of the self-energy with re- 
spect to frequency is momentum independent (in contrast 
to the static part), yielding a k independent QP residue. 
Hence, our approach goes beyond these methods, at least 
in conjunction with the cluster-DMFT approach. We 
hope that having an economic impurity solver based on 
SBs will allow us to study larger cluster sizes than feasible 
with other methods, and most importantly help us un- 
derstanding the low-energy physics emerging from these 
cluster dynamical mean- field theories. 

Finally, we limited our study to slave bosons which 
do not mix the particle number. The extension to 
full charge-rotational invariance and superconductivity 
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is possible (sec Refs. [14,22] in the single-orbital case), 
and will be presented in a separate paper. In this con- 
text, the slave boson method will incorporate the SU{2) 
charge symmetry and its extension away from half fill- 
ing considered by Wen and Lee^^ and the rotationally- 
invariant slave-boson formalism can serve as a power- 
ful tool for interpreting the low-energy physics emerging 
from plaquette-CDMFT studies of this issue. 
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APPENDIX A: SINGLE-ORBITAL CASE AND 
CONNECTION WITH PREVIOUS WORK 

Here, we briefly consider the single-orbital case 
(M=2), which also allows to make contact with Refs. [13, 
14] . These authors introduced in this case a rotationally- 
invariant formalism, with the calculation of response 
functions associated with the saddle point as their main 
motivation. For M=2, the following local basis set can 
be considered (whether or not H\oc is diagonal in this 
basis): 



= 
A^ = 1 
N = 2 



|0), 

[£)) = d|d/|0) 



(Al) 



Hence, wc introduce the following bosons (not mixing 
sectors with different particle numbers, i.e., not consid- 
ering superconducting states): 



<i>00 = <I>E 



0Ti = 



(A2) 



Up to normalizations, the bosons p^^, introduced in 
Ref. 13 correspond to 0^^/- In contrast, the stan- 
dard Kotliar-Ruckcnstcin^^ formalism introduces only 
two bosons p1 in the one-particle sector. The represen- 



tatives (27) of the physical states read here: 
|0) = 4Kac) 



1^) = 4E'^L'/]'|vac) 



\D) = (p\jdj dl\va.c) 
and the constraints (28,29) read: 



(A3) 



(A4) 

(A5) 
(A6) 
(A7) 



Not including, for simplicity, the square-root normaliza- 
tions in (37), needed however in order to insure a correct 
U =0 limit at saddle point, the 'simplest' expression (32) 
of the electron creation operators read: 



1 




aa' 


f 




"E^^a'/'aa 
cr 


f\fi 


(7 








Ki ■ 



4 = 7fE['^U^-(-iM'^d4 • 



(A8) 
(A9) 



Apart from the motivations of Ref. [13] (associated with 
fluctuations and response functions), the usefulness of the 
rotationally-invariant scheme in the single-orbital case 
can be demonstrated on a toy model consisting of a one- 
band Hubbard model with a magnetic field, purposely 
written in the direction (i.e.. in the form hd}^di +h.c, 
analogous to a hybridization) . Although the direction of 
the field should not matter, a direct application of the 
standard Kotliar-Ruckenstein formalism is impossible in 
that case. The rotationally-invariant formalism can be 
shown to lead to the correct saddle point, independently 
of the spin-quantization axis. 



APPENDIX B: DERIVATION OF EQ. (29) 

In this section, wc show that the physical states of the 
form (27) are exactly those selected by the constraints 
(29) and (28). First, it is easy to check that states of the 
form (27) do satisfy these constraints. Indeed, let us act 



on the state |C) 
(29). The l.h.s leads to: 

1 



''Cm 



vac) ®\m)f and with 



flU'\c) = 




Cmlvac) (8)/J,/a'|TO)/ 

m\fl fa'\m) (/-l^^jvac) O |m')/ 



(Bl) 
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When acting with the r.h.s, only the term A=C and n=m 
gives a non-vanishing contribution, hence: 



A nn' 

-/^IZH/l/a'l"'>4n'|vac)(8)|n)/ . 

^ nn' 

We now prove that (29) are sw/Jjcient conditions, which 
is a bit more difficult. Since (28) excludes states with 
more than one boson, it is enough to consider a general 
state of the form: 

|C; W)=Y, Wp<i '/'Cplvac) ® \q)f . (B2) 

pq 

and to show that (29) implies Wpq ex 5pq. Acting on this 
state with each term in the constraint (29) yields for the 
l.h.s: 

fi fa'\C; W) =Y^ Wpq 4p|vac) ® flUq)f 

pq 

= J2 4p|vac) \r)f J2 WpqMUq) .(B3) 

pr q 

Let US now act with the r.h.s. Only the terms with A=C 

and n=p contribute, leading to: 



EE-^Jin'-^An (n|/t/„,|n')|C;T^) 

A nn' 

= E Wpq 4n'kac) ® \q)f{p\fiU'\n') 

pqn' 

= E-^Cplvac) ® \r)fY.WAq\fifc.'\p) ,(B4) 

pr q 

where the last expression comes from a change of indices 
n' ^p,q-fr,p-f q. 

We see that the constraint is satisfied provided that 
the following identity holds, for all orbital indices aa' 
and all states p, r: 

J2Wp,{r\fifa'\q)=J2W,r{q\fya'\p) . (B5) 



Let us first look at the case a^a' , which reads: 

ra Wpr = Pa Wpr ■ 

Hence Wpr=0 unless Pa=ra for all a, so that: 



pq 



^pq 



Substituting this into (B5), we obtain: 

Wp {r\flfa' \p) = Wr {rlflfa' \p) ■ 



(B6) 
(B7) 
(B8) 



Thus, if r and p are related by a move of a QP from 
one state to another (a transposition of two occupation 

numbers), then Wp=Wr- Moreover, two Fock states in 
the same sector Hn of the Hilbert space are related by 
a permutation of the occupied states, which can be de- 
composed in a product of transpositions. Hence, Wp is a 
constant for p e Hpf, and WpqfxSpq as claimed. 



APPENDIX C: PHYSICAL CREATION 
OPERATOR 

1. Proximate expression 

First let us note that there is a systematic route to 
find the expression for d, which consists in writing the 
operator as : 



di 



{A\dl\B) \A}{B\ 



AB 

E 

AB neHA,meHB 



J2 {A\dl\B)cl,l^4>BmXi^, (CI) 



with X/^j^~\n) f{m\f in usual Hubbard notations. X^^^ is 
obviously noi just a one-particle operator f^, even when 
restricted to the sectors of interest in the above formula, 
since the states n and m can differ in many places. How- 
ever, because any transposition of two QPs, when acting 
on a physical state, can be replaced by a corresponding 
operation on bosons using the constraint (29), and be- 
cause any product of bosonic operators which cannot be 
reduced to a quadratic form will produce a state which is 
out of the physical subspace, the physical operator must 
in the end take the form: 



^l = EEE^. 

13 AB nm 



t 



(C2) 



One can solve for the coefficients C^^(q:,/J), requesting 
proper action on the physical states. 

In this section however, we restrict ourselves to proving 
that Eq. (32) does the job, i.e., that 



{A\di\B){n\fl\m) 



4= y , 

,,knn. ^N^{M-Ns) 

satisfies 



il\B)=Y.{A\dl\B) \A) 



We start by proving the formula : 

E E {n\fl\p)fl\p) = {N + l)\n) 
peHN 



(C3) 



(C4) 



(C5) 



where the sum over p runs over the basis of the subspace 
Hn (states with QPs) of the Fock space. First, we 
have in general : 



EE 

/3 peHN 



{n\f;\p) fl\p) 



n'eHN+i 



(C6) 



but, because is the creation operator (it connects one 
basis state to only one another) : 



an' 



E {n\4\p}{n'\fl\p)c<6nn' E 

(C7) 
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Wc now use the fact that the tensor is invariant 
(it has the same expression in every basis) and use 
the notations introduced for Eq. (40): U is a uni- 
tary transformation of the one QP states and U is 
the corresponding transformation in the Fock states. 
{n\f^^\m)=U0i3'l(*n,{n'\f^^,\m')l(mm' and we have : 

^ ^ I ,2 

= Y Y Y U0i3'U00M*„,U„n„ 



l3,meHtf 0'm'€Hff 0"m"€Hff 
n'SHjv + l i"eHN + l 



^ ^mm'^mm" 



(n'|4'|m')(n"|/^„|m") 



"I ft 



EE E Kn'i^nn"{n'\fl,\p){n'Vh\Pr 

0' peHp, n',n"eHN+i 

2 



J2 |((Wn)|4|m) 



(C8) 



Moreover, any couple of c;lements of the basis of the Fock 
state can be connected by a W transformation (with a U 
that permutes the one QP basis state), therefore o„=a is 
a independent of n. a can then be determined by sum- 
ming (C5) over n : 



E 



E E \(^\f>) 



f3 "e«N + i 



E hE4//'I") 
E (^+1) ' 



(C9) 



leading to a=N + 1. This completes the proof of (C5). 

It is now simple to compute the action of (32) : acting 
on \C)=^^^^ J2p 0Cp|vac) ig) \p) / with this operator, only 
the term B=C,m=p contributes, and we get : 

4ic) = ^ 

^/Dc{Nc + 1)(M - Nc) 

13 P&Hnc 

Y,{A\dl\C)<l,\jvac)®\n) 
(CIO) 



Nc + l 



Dc{M - Nc) 



= J2i^\dl\C)\C) , 

A 

which is identical to (31). 

2. Improved expression 



In this section we present arguments for the improved 
formula used in this paper. First, it is useful to define the 



"natural orbitals" (NO) basis as the basis which diago- 
nalizes the quasiparticle and quasihole density matrices 
corresponding to the average constraint, which is given 

by 



E '^in^Am(m|/^//3|n) (Clla) 



Anm 



^ E ^*An'pAm{m\Ufl\n) (CI lb) 



Anm 



An 

Let us denote by ^Aj |A) the eigenvalues and eigenvectors 
of those matrices : 

^<^:l=Y^x{a\X){X\/3) , \X)=Y{a\X)\a) , (C12) 

A A 

which is equivalent to use the NO quasiparticle operator 
'ip\ such that : 



(C13) 



To be fully explicit, we can consider the particular basis 
transformation (in the notations of the section above) 
fa—^ctxi^x which rotates to the NOs, and the corre- 
sponding rotation on the bosons: <pAn=^{U)nn'^An' ■ 
The rotation matrix is: 

f/aA = (a|A) , (C14) 

and in the NO basis: 

Y n*^nnBm{m\tjjltjjij,\n) = 5a/h E ^^"^^"'^^ 



Anm 



An 



and: 



.(p)r 



= Sx^^xii^^An}) , (C15) 



J2^aX^x[U^]x0 . (C16) 



The idea is to generalize the Kotliar-Ruckenstein nor- 
malization factor in the NO basis, where the QP density 

being diagonal, its probabilistic interpretation is more 
transparent. Hence, the improved expression of d reads: 



~" A Is'nm VU{^An}){l " Ui^An})) 

(C17) 

Note that the formal square-root normalisation, i.e., 
1/ ^Na{M — Nb), does not appear in this representa- 
tion. We can now rotate back to the generic basis we 
started from and use the gauge invariance, leading to : 

rAn/^ mj,t J, \^ (/3|A)(A|7) ^ 
dl = 2^ Csmia, /?) (p'An't'Bm l^—rf ,._,J l 
AB,nm,0j A VsAU t.X) 

= C^;j,(a,/3)^i„</.5„(/3|[A(^')AW]-^|7)/t ,(C18) 
with C^:^{a,P) = {A\di\B}(n\fl\m) .(C19) 



19 



Hence this yields the following form for the i?-matrix: 



grange multipliers: 



(C20) 

In the actual implementation of the saddle point calcu- 
lations, the explicit use of both the quasiparticle and the 
quasihole density matrices has been utilized (i.e., not us- 
ing their relation) . Although at convergence the different 
representations yield the same values, writing the equa- 
tions via both, particle and hole density matrix, appears 
to be necessary within the minimization cycle. This is 
due to the fact that the derivatives with respect to the 
slave bosons have to be symmetric, however when using 
only the particle density matrix (or its eigcnsystem de- 
composition) for instance the derivative with respect to 
the empty boson vanish, although this one exists in the 
KR case. In the end an even more symmetrized form, 
i.e., i(A(p)A('') + AC'^A^p)) was used for the square root 
in eq. (C18). Thus defining the following matrix; 



P) , (C21) 



the electron operators are written as 



AB nm ■y/3 

= E</34 (C22) 

AB nm 7/3 

= J2Ra0f0 , (C23) 



and correspondingly independently written the elements 
of the i?,i?'^-matrices read 



i?M„/3= E CBmic^^l)<^^Bm^AnM0^ (C24) 

R^[<PU = = E ^BliPn) <l>\MmM^a. {C25) 

AB,nm,j 



APPENDIX D: DETAILS ON THE 
SADDLE-POINT EQUATIONS AND THEIR 
NUMERICAL SOLUTION 



The saddle-point equations for T=0 are obtained by 

performing the partial derivatives with respect to all the 
variables, i.e., condensed slave-boson amplitudes and La- 



An 



dX 
dfl 



(Dl) 



= {flM- E 'pL' i^ldidfiW) ^An (D2) 



A,nn' 



-EA«/3E("^l'^a^/3|"'>^Cn' (D3) 

a/3 n' 

E /kj + E ^CB^Bm + ^0<Pcm 

■■'Cm kj "^Cm B 

- E E ("I 1 "fiCn (D4) 

a/3 n 

With (/^//3)=E/k.(«kk,>K,|/3) . (D5) 

The Ekj are the eigenvalues (with band index j) of the QP 
matrix (R^((p)£(k)R((^) + A) with corresponding eigen- 
vector |fkj), while /kj denotes the occupation number of 
the state |fkj) for a given total number of particles, to 
be evaluated by standard fc-integration techniques (e.g. 
tetrahedron method, Gaussian smearing, etc.). 



1. Some slave-boson derivatives 

a. Eigenvalues and R matrices. The derivatives of 
the eigenvalues with respect to the slave bosons, i.e., 
may be performed pertubatively: 



de 



kj 



dip 



= (t'kj 



d_ 
dip 
5Rt 



(Rt£(k)R + A) 



J^ki 
f)Ti 

— s(k)R + Rts(k)- 



^3 



dp 



+ ER^-(fc)l«)^(/3 



a/3 



E 

al3 



('^k,-|a)^(/?|e(k)R|^k,-) 
+K,-|R^£(k)|a)^(/3K,) 



.(D6) 



The therefore needed explicit expressions for the deriva- 
tives of the i?,i?^-matrices read as follows (using 
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eqs.(C24,C25)): 



2. Mixing 



dR 



d^Cm 

dRafi 

d^pT 



AB ,nn' 



Cm AB ,nn' ,j 

^analogous for [^^]„/3^ 

b. The M matrix. As it is seen, the derivatives in- 
volve the derivative of the M matrix (C21). This deriva- 
tive is computed as follows. Lets first write M as 



7/3 



(7|K-V2|/J) _ 



(D7) 



What we are looking for is the derivative of K~^/^ with 
respect to the SBs. In order to get access to this quantity 

we use the identity 



(S^K-V2)K-i/2 + K-i/2(a^K-i/2) 



a^K-i (D8) 
Y , (D9) 



with X=9^K-i/2 and Y=9^K-i. We then apply P 
which transforms K to its eigensystem. This yields 



X'L + LX' = Y' 



(DIG) 



where the prime denotes that the quantities de- 
fined above are expressed in that eigensystem, and 
L=P'I'K~^/^P. Since in the eigensystem K~^/^ is di- 
agonal, i.e., L is, the last equation can be written in 
components and X' determined: 



A', 



(Dll) 



Backtransforming to X=PX'P'I'=9^K~^/^ yields the 
desired derivative of K and subsequently of M. 
To perform the described computation we need to 
know d:pK~^ in eq. (D8), however this quan- 
tity may be straightforwardly calculated when start- 
ing from the identity KK~^^^K~-'^/^=l, resulting in 



-K-i(5yK)K- 



In order to solve the saddle-point equations, a method 
to deal with a system of nonlinear equations F as a func- 
tion of the variables (slave bosons, Lagrange multipliers) 
X has to be utilized: 



F(x) = 



(D12) 



In the present work we tested several quasi-Newton tech- 
niques (e.g. Broyden"'''^, modified Broyden"'''*, etc.) to 
handle this numerically. Thereby from a starting guess 
for X the variables are updated via 



X 



since we want F 



(m+l) 



^^M+Jrm)F(m) , (D13) 

to be zero to linear order. The 



jacobian J is here defined as follows 



(D14) 



and is not calculated exactly (this would involve sec- 
ond derivatives and would lead to the Newton-Raphson 
method) but is computed at each step m via formulae 
which dictate several constraints on how J should evolve. 
In our numerical implementation we found the modified 
Broyden scheme to be well suited for the so far investi- 
gated applications Note that there is usually no need for 
explicitly fixing the gauge for the numerical solution of 
the saddle-point equations. The initial amplitudes of the 
variational parameters, i.e., slave bosons and lagrange 
multipliers, at the start of the iteration, together with 
the choice of the atomic basis \A), always ensured proper 
convergence to one of the family of solutions within our 
implementation. 



* Electronic address: Frank. Lechermann@physnet . 
uni-hamburg . de 

^ G. Kotliar, in Strongly Interacting Fermions and High-Tc 
superconductivity, edited by B. Doucot and J. Zinn-Justin 
(Elsevier, 1995), Les Houches, Session LVI, p. 197. 

^ M. Imada, A. Fujimori, and Y. Tokura, Rev. Mod. Phys. 
70, 1039 (1998). 

^ A. Damascelli, Z. Hussain, and Z.-X. Shen, Rev. Mod. 
Phys. 75, 473 (2003), cond-mat/0208504. 

* M. R. Norman et al, Nature 392, 157 (1998). 

® M. Le Tacon, A. Sacuto, A. Georges, G. Kotliar, Y. Gallais, 



D. Colson, and A. Forget, Nature Physics 2, 537 (2006), 
arXiv:cond-mat/0603392. 
® T. Maier, M. Jarrell, T. Pruschke, and M. H. Hettler, Rev. 
Mod. Phys. 77, 1027 (2005). 

G. Kotliar, S. Y. Savrasov, K. Haule, V. S. Oudovenko, 
O. ParcoUet, and C. A. Marianetti, Rev. Mod. Phys. 78, 

865 (2006). 

^ A. M. S. Trcmblay, B. Kyung, and D. Scncchal, ,]. Low 

Temp. Physics 32, 424 (2006). 
^ A. Georges, G. Kotliar, W. Krauth, and M. J. Rozenberg, 

Reviews of Modern Physics 68, 13 (1996). 



21 



^° D. Senechal and A.-M. S. Tremblay, Phys. Rev. Lett. 92, 

126401 (2004). 

M. CivcUi, M. Caponc, S. S. Kancharla, O. Parcollet, and 
G. Kotliar, Phys. Rev. Lett. 95, 106402 (2005). 
O. Parcollet, G. Biroli, and G. Kotliar, Phys. Rev. Lett. 
92, 226402 (2004). 
13 T. Li, P. Wolfle, and P. J. Hirschfeld, Phys. Rev. B. 40, 
6817 (1989). 

" R. Fresard and P. Wolfle, Int. J. Mod. Phys. B 6, 685 
(1992). 

C. Attaccalite and M. Fabrizio, Phys. Rev. B 68, 155117 
(2003). 

1^ M. Ferrero, Ph.D. thesis, SISSA-Trieste (2006). 
^'^ G. Kotliar and A. E. Ruckenstein, Phys. Rev. Lett. 57, 
1362 (1986). 

1* R. Fresard and G. Kotliar, Phys. Rev. B 56, 12909 (1997). 
1® J. Biinemann, W. Weber, and F. Gebhard, Phys. Rev. B 

67, 6896 (1998). 
^° B. R. Trees, A. J. Fedro, and M. R. Norman, Phys. Rev. 

B 51, 6167 (1995). 

21 X. Dai, G. Kotliar, and Z. Fang, cond-mat/0611075 (2006). 

22 B. R. Bulka and S. Robaszkiewicz, Phys. Rev.B 54, 13138 
(1996). 

23 J. H. Shim, K. Haule, and G. Kothar, Nature 446, 513 
(2007). 

2* Y. Ono, M. Potthoff, and R. Bulla, Phys. Rev. B 67, 

035119 (2003). 
25 M. Rozcnbcrg, Phys. Rev. B 55, R4855 (1997). 

T. Pruschke and R. Bulla, Eur. Phys. J. B 44, 217 (2005). 
2'' K. Inaba and A. Koga, Phys. Rev. B 73, 155106 (2006). 
2* J. E. Han, M. Jarrefl, and D. L. Cox, Phys. Rev. B 58, 

R4199 (1999). 

2^^ Y. Song and L.-J. Zou, Phys. Rev. B 72, 085114 (2005). 
3° A. Koga, N. Kawakami, T. M. Rice, and M. Sigrist, Phys. 
Rev. Lett. 92, 216402 (2004). 

31 A. Riiegg, M. Indergaard, S. Pilgram, and M. Sigrist, Eur. 
Phys. J. B 48, 55 (2005). 

32 L. de' Medici, A. Georges, and S. Biermann, Phys. Rev. B 
72 (2005), cond-mat/0503764. 

33 V. I. Anisimov, I. Nekrasov, D. Kondakov, T. M. Rice, and 



M. Sigrist, Eur. Phys. J. B 25, 191 (2002). 
3* A. Licbsch, Phys. Rev. B 70, 165103 (2004). 
35 M. Ferrero, F. Bccca, M. Fabrizio, and M. Capone, Phys. 

Rev. B 72, 205126 (2005). 
3^ H. Monien, N. Elstner, and A. J. Millis, cond-mat/9707051 

(1997). 

3^^ G. Moeller, V. Dobrosavljevic, and A. E. Ruckenstein, 

Phys. Rev. B 59, 6846 (1999). 
3® A. Fuhrmann, D. Heilmann, and H. Monien, Phys. Rev. B 

73, 245118 (2006). 
3" S. S. Kancharla and S. Okamoto, cond-mat/0703728 

(2007). 

'"^ M. Capone, private communication (2007). 
"1 M. Ferrero, P. S. Cornagha, and et al, unpublished (2007). 
■'2 G. Biroli, O. Parcollet, and G. Kotliar, Phys. Rev. B 69, 
205108 (2004). 

"3 G. Kotliar, S. Y. Savrasov, G. Palson, and G. BiroU, Phys. 

Rev. Lett. 87, 186401 (2001). 

J. E. Hirsch, Phys. Rev. B 31, 4403 (1985). 
*5 N. Furukawa and M. Imada, J. Phys. Soc. Jpn. 61, 3331 

(1992). 

G. Fano, F. Ortolani, and A. Parola, Phys. Rev. B 42, 
6877 (1990). 

E. Dagotto, A. Moreo, F. Ortolani, D. Poilblanc, and J. Ri- 
era, Phys. Rev. B 45, 10741 (1992). 

T. Kashima and M. Imada, J. Phys. Soc. Jpn. 70, 2287 
(2001). 

D. Rohe and W. Mctznor, Phys. Rev. B 71, 115116 (2005). 
5° T. Maier, M. Jarrell, T. Pruschke, and M. H. Hettler, Rev. 
Mod. Phys. 77, 1027 (2005). 

51 X.-G. Wen and P. A. Lee, Phys. Rev. Lett. 76, 503 (1996). 

52 M. Fabrizio, unpubhshod (2007). 

53 C. G. Broydon, Math. Comp. 19, 577 (1965). 

5* D. Vandcrbilt and S. Louie, Phys. Rev. B 30, 6118 (1984). 
55 In this paper, the chemical potential is usually included in 

the one-body part of the local hamiltonian. 
5® Note however that both sides of (24) have identical matrix 

elements between physical states 



